library(nnet)
library(rms)
data<-read.csv(file.choose())


# 3. 处理类别变量（自动转换哑变量）
data$category<-as.factor(data$category)
data$category<-relevel(data$category, ref = "3")
data$education<-as.factor(data$education)
data$education<-relevel(data$education, ref = "4")
data$hypertension<-as.factor(data$hypertension)
data$hypertension<-relevel(data$hypertension, ref = "1")
data$storke<-as.factor(data$storke)
data$storke<-relevel(data$storke, ref = "1")
data$Memory.disorder<-as.factor(data$Memory.disorder)
data$Memory.disorder<-relevel(data$Memory.disorder, ref = "1")
data$three.high<-as.factor(data$three.high)
data$three.high<-relevel(data$three.high, ref = "1")
data$drink<-as.factor(data$drink)
data$drink<-relevel(data$drink, ref = "1")
data$BMI<-as.factor(data$BMI)
data$BMI<-relevel(data$BMI, ref = "4")
data$depression<-as.factor(data$depression)
data$depression<-relevel(data$depression, ref = "3")
data$sleep<-as.factor(data$sleep)
data$sleep<-relevel(data$sleep, ref = "3")
data$disability<-as.factor(data$disability)
data$disability<-relevel(data$disability, ref = "1")
data$hearing<-as.factor(data$hearing)
data$hearing<-relevel(data$hearing, ref = "5")




# 重新创建 datadist
dd <- datadist(data)
options(datadist = "dd")

data1<-subset(data,category %in% c(1,3))
data1$category <- ifelse(data1$category == 1, 1, 0)  # 重新编码为 0 和 1
data2<-subset(data,category %in% c(2,3))
data2$category <- ifelse(data2$category == 2, 1, 0)  # 重新编码为 0 和 1
fit1 <- lrm(category ~  bl_ua+education+hypertension+hearing+Memory.disorder+BMI+
             three.high+depression+sleep+disability+storke+drink, data = data1)
fit1$coefficients<-c(7.615,0.687,-1.219,-0.477,0.980,
                     -0.573,-4.086,-0.418,-0.477,3.722,
                     -2.160,-2.794,-0.558,0.696,
                     -2.130,-5.063,-3.269,-4.280,
                     2.099,1.847,-1.700,1.363,
                     -3.300)

nom <- nomogram(fit1, fun = plogis, funlabel = paste("Probability of Class 1"),fun.at = seq(0.1, 0.9, by=0.8),lp.at = seq(-10, 70, by=5))
plot(nom)

fit2 <- lrm(category ~  bl_ua+education+hypertension+hearing+Memory.disorder+BMI+
              three.high+depression+sleep+disability+storke+drink, data = data2)
fit2$coefficients<-c(1.750,0.494,-4.466,-2.353,-0.563,-0.991,
                     -1.005,2.202,1.786,4.695,0.313,-3.599,
                     -0.850,0.074,-0.502,-2.311,-1.551,-2.028,
                     2.673,3.051,-0.080,0.259,-0.748)
nom1 <- nomogram(fit2, fun = plogis, funlabel = paste("Probability of Class 2"),fun.at = seq(0.1, 0.9, by=0.2))
plot(nom1)





####随机森林
library(rpart.plot)
library(rpart)
library(randomForest)
library(xgboost)
library(dplyr)
library(caret)
library(pROC)
library(skimr)
library(DataExplorer)
library(ggplot2)
library(reshape2)
library(dplyr)
library(ggpubr)
library(showtext)
library(partykit)

##拆分数据
set.seed(1234)#固定划分数据结果
trainIndex <- createDataPartition(data$category, p = 0.7, list = FALSE)
train <- data[trainIndex, ]
test <- data[-trainIndex, ]
##构建模型
set.seed(12345)#固定交叉验证结果，保证结果的可重复性
rf_model <- randomForest(x=train[,c(2:13)],
                         y=train$category,
                         data = train, 
                         method="class",
                         ntree = 1000, #决策树颗数
                         importance = T #输出变量重要性
)
##变量重要性
importance(rf_model)
# 设置图形布局为 1 行 2 列
par(mfrow = c(1, 2))  # 1 行 2 列
# 绘制 Mean Decrease Accuracy (MDA)
varImpPlot(rf_model,  type = 1)
# 绘制 Mean Decrease Gini (MDG)
varImpPlot(rf_model,  type = 2)
# 恢复默认图形布局
par(mfrow = c(1, 1))  # 恢复为 1 行 1 列
#################
## 预测

## 训练集预测概率
train_pred <- predict(rf_model, newdata = train, type = "prob")

## 训练集预测分类
trainpredlab <- predict(rf_model, newdata = train, type = "class")

## 训练集混淆矩阵
confusionMatrix(data = trainpredlab,reference = train$category, mode = "everything")

## 测试集预测概率
test_pred <- predict(rf_model, newdata = test, type = "prob")

## 测试集预测分类
testpredlab <- predict(rf_model, newdata = test, type = "class")

## 测试集混淆矩阵
confusionMatrix(data = testpredlab,reference = test$category, mode = "everything")






